home *** CD-ROM | disk | FTP | other *** search
/ AOL File Library: 2,801 to 2,900 / aol-file-protocol-4400-2801-to-2900.zip / AOLDLs / C++ Files Library / Graphic Gems I, II & III (C_C++) / Graphics Gems C Code.sea / GemsII / Peano / peano.c < prev    next >
Text File  |  1992-06-16  |  4KB  |  232 lines

  1. /* peano.c */
  2.  
  3. /* copyright Ken Musgrave */
  4. /* June 1986 */
  5.  
  6. /*
  7.  * space-filling peano curve algorithm.
  8.  * fills n-space with a 1-D peano curve.
  9.  * 
  10.  * the algorithm is utterly incomprehensible,
  11.  * so expect a paucity of sensible comments in this code.
  12.  * it consists largely of bit-wise logical operations,
  13.  * and is therefore quite inpenetrable.
  14.  * but it works.
  15.  */
  16.  
  17.  
  18. #include "types.h"
  19.  
  20.  
  21.  
  22. /*
  23.  * determine the n-space coordinate of "point" on the peano curve
  24.  */
  25. peano(coord, point)
  26. vector          coord;
  27. int             point;
  28. {
  29.     int             i;
  30.  
  31.     zero(sigma);        /* initialize necessary arrays */
  32.     zero(tilde_sigma);
  33.     zero(tilde_tau);
  34.  
  35.     build_rho(point);
  36.     for (i=0; i<precision; i++)
  37.         J[i] = principal_pos(rho[i]);
  38.     build_sigma();
  39.     build_tau();
  40.     build_tilde_sigma();
  41.     build_tilde_tau();
  42.     build_omega();
  43.     build_alpha();
  44.  
  45.     v_convert(alpha, coord);
  46.  
  47. } /* peano() */
  48.  
  49.  
  50. /*
  51.  * build "rho" array
  52.  */
  53. build_rho(point)
  54. int             point;
  55. {
  56.     int    i, mask=bytemask;
  57.  
  58.     for (i=0; i<precision; i++) {
  59.         rho[precision - i - 1] = (point & mask) >> (i * dimensions);
  60.         mask <<= dimensions;
  61.     }
  62.  
  63. } /* build_rho() */
  64.  
  65.  
  66. /*
  67.  * find principal position of "a_byte" 
  68.  */
  69. principal_pos(a_byte)
  70. byte            a_byte;
  71. {
  72.     int             nth_bit, i=1;
  73.  
  74.     nth_bit = a_byte & 0x01;
  75.     for (i=1; i<dimensions; i++) {
  76.         if (((a_byte & bitmask[dimensions - i - 1]) >> i) != nth_bit)
  77.             return (dimensions - i);
  78.     }
  79.  
  80.     return (dimensions);    /* all bits are the same */
  81.  
  82. } /* principal_pos() */
  83.  
  84.  
  85. /*
  86.  * build "sigma" array
  87.  */
  88. build_sigma()
  89. {
  90.     int             i, bit;
  91.  
  92.     for (i=0; i<precision; i++) {
  93.         sigma[i] |= rho[i] & bitmask[0];
  94.         for (bit=1; bit<dimensions; bit++) {
  95.             sigma[i] |= (rho[i] & bitmask[bit])
  96.                 ^ ((rho[i] & bitmask[bit - 1]) >> 1);
  97.         }
  98.     }
  99.  
  100. } /* build_sigma() */
  101.  
  102.  
  103. /*
  104.  * build "tau" array 
  105.  */
  106. build_tau()
  107. {
  108.     int             parity, bit, i, j;
  109.     byte            temp_byte;
  110.  
  111.     for (i=0; i<precision; i++) {
  112.         parity = 0;
  113.                         /* complement nth bit */
  114.         if (sigma[i] & bitmask[dimensions - 1])
  115.             tau[i] = sigma[i] - 1;    /* nth bit was 1 */
  116.         else
  117.             tau[i] = sigma[i] + 1;    /* nth bit was 0 */
  118.  
  119.         for (j=0; j<dimensions; j++)    /* find parity */
  120.             if (tau[i] & bitmask[j])
  121.                 parity++;
  122.  
  123.         if (ODD(parity)) {    /* complement principal bit */
  124.             bit = J[i] - 1;    /* get index of principal bit */
  125.                     /* get bit in question */
  126.             temp_byte = tau[i] & bitmask[bit];
  127.             tau[i] |= bitmask[bit];    /* set the bit to 1 */
  128.             tau[i] &= ~temp_byte;    /* assign complement */
  129.         }
  130.     }
  131.  
  132. } /* build_tau() */
  133.  
  134.  
  135. /*
  136.  * build "tilde_sigma" array 
  137.  */
  138. build_tilde_sigma()
  139. {
  140.     int             i, shift=0;
  141.  
  142.     tilde_sigma[0] = sigma[0];
  143.     for (i=1; i<precision; i++) {
  144.         shift += J[i - 1] - 1;
  145.         shift %= dimensions;
  146.         tilde_sigma[i]= RT_CSHFT(sigma[i], shift, dimensions, bytemask);
  147.     }
  148.  
  149. } /* build_tilde_sigma() */
  150.  
  151.  
  152. /*
  153.  * build "tilde_tau" array 
  154.  */
  155. build_tilde_tau()
  156. {
  157.     int             i, shift=0;
  158.  
  159.     tilde_tau[0] = tau[0];
  160.     for (i=1; i<precision; i++) {
  161.         shift += J[i - 1] - 1;
  162.         shift %= dimensions;
  163.         tilde_tau[i] = RT_CSHFT(tau[i], shift, dimensions, bytemask);
  164.     }
  165.  
  166. } /* build_tilde_tau() */
  167.  
  168.  
  169. /*
  170.  * build "omega" array 
  171.  */
  172. build_omega()
  173. {
  174.     int             i;
  175.  
  176.     omega[0] = 0;
  177.     for (i=1; i<precision; i++)
  178.         omega[i] = omega[i - 1] ^ tilde_tau[i - 1];
  179.  
  180. } /* build_omega() */
  181.  
  182.  
  183. /*
  184.  * build "alpha" array 
  185.  */
  186. build_alpha()
  187. {
  188.     int             i;
  189.  
  190.     for (i=0; i<precision; i++)
  191.         alpha[i] = omega[i] ^ tilde_sigma[i];
  192.  
  193. } /* build_alpha() */
  194.  
  195.  
  196. /*
  197.  * initialize "array" to zeros 
  198.  */
  199. zero(array)
  200.  
  201.     r_array         array;
  202.  
  203. {
  204.     int             i;
  205.  
  206.     for (i=0; i<precision; i++)
  207.         array[i] = 0;
  208.  
  209. } /* zero() */
  210.  
  211.  
  212. /*
  213.  * convert "alpha" array into n_space coordinate vector 
  214.  */
  215. v_convert(alpha, coord)
  216. r_array         alpha;
  217. vector          coord;
  218. {
  219.     int             i, j, bit, a_bitmask=1;
  220.  
  221.     for (i=0; i<dimensions; i++) {
  222.         coord[i] = 0;
  223.         bit = precision;
  224.         for (j=0; j<precision; j++)    /* extract each bit of coord
  225.                          * i */
  226.             coord[i] |= ((alpha[j] & a_bitmask) << --bit) >> i;
  227.         a_bitmask <<= 1;
  228.     }
  229.  
  230. } /* v_convert() */
  231.  
  232.